import sys
from datetime import datetime

from astropy.io import fits as pyfits
from astropy.visualization.interval import ZScaleInterval
import matplotlib.pyplot as plt
from pylab import *
import numpy as np
import PIL.Image as pil
from astropy.wcs import WCS

#https://ps1images.stsci.edu/ps1image.html
#from __future__ import print_function
import numpy
from astropy.table import Table
import requests
from PIL import Image
from io import BytesIO


def xprint_time():

    from datetime import datetime
    # 获取当前时间
    now = datetime.now()
    # 打印当前时间
    print("Current time:", now.strftime("%Y-%m-%d %H:%M:%S"))

def getimages(ra, dec, size=500, filters="grizy"):
    """Query ps1filenames.py service to get a list of images

    ra, dec = position in degrees
    size = image size in pixels (0.25 arcsec/pixel)
    filters = string with filters to include
    Returns a table with the results
    """

    service = "https://ps1images.stsci.edu/cgi-bin/ps1filenames.py"
    url = ("{service}?ra={ra}&dec={dec}&size={size}&format=fits"
           "&filters={filters}").format(**locals())
    table = Table.read(url, format='ascii')
    return table


def geturl(ra, dec, size=500, output_size=None, filters="grizy", format="jpg", color=False):
    """Get URL for images in the table

    ra, dec = position in degrees
    size = extracted image size in pixels (0.25 arcsec/pixel)
    output_size = output (display) image size in pixels (default = size).
                  output_size has no effect for fits format images.
    filters = string with filters to include
    format = data format (options are "jpg", "png" or "fits")
    color = if True, creates a color image (only for jpg or png format).
            Default is return a list of URLs for single-filter grayscale images.
    Returns a string with the URL
    """

    if color and format == "fits":
        raise ValueError("color images are available only for jpg or png formats")
    if format not in ("jpg", "png", "fits"):
        raise ValueError("format must be one of jpg, png, fits")
    table = getimages(ra, dec, size=size, filters=filters)
    url = ("https://ps1images.stsci.edu/cgi-bin/fitscut.cgi?"
           "ra={ra}&dec={dec}&size={size}&format={format}").format(**locals())
    if output_size:
        url = url + "&output_size={}".format(output_size)
    # sort filters from red to blue
    flist = ["yzirg".find(x) for x in table['filter']]
    table = table[numpy.argsort(flist)]
    if color:
        if len(table) > 3:
            # pick 3 filters
            table = table[[0, len(table) // 2, len(table) - 1]]
        for i, param in enumerate(["red", "green", "blue"]):
            url = url + "&{}={}".format(param, table['filename'][i])
    else:
        urlbase = url + "&red="
        url = []
        for filename in table['filename']:
            url.append(urlbase + filename)
    #print(url)
    return url


def getcolorim(ra, dec, size=100, output_size=None, filters="grizy", format="jpg"):
    """Get color image at a sky position

    ra, dec = position in degrees
    size = extracted image size in pixels (0.25 arcsec/pixel)
    output_size = output (display) image size in pixels (default = size).
                  output_size has no effect for fits format images.
    filters = string with filters to include
    format = data format (options are "jpg", "png")
    Returns the image
    """

    if format not in ("jpg", "png"):
        raise ValueError("format must be jpg or png")
    url = geturl(ra, dec, size=size, filters=filters, output_size=output_size, format=format, color=True)
    r = requests.get(url)
    im = Image.open(BytesIO(r.content))
    return im


def getgrayim(ra, dec, size=100, output_size=None, filter="g", format="jpg"):
    """Get grayscale image at a sky position

    ra, dec = position in degrees
    size = extracted image size in pixels (0.25 arcsec/pixel)
    output_size = output (display) image size in pixels (default = size).
                  output_size has no effect for fits format images.
    filter = string with filter to extract (one of grizy)
    format = data format (options are "jpg", "png")
    Returns the image
    """

    if format not in ("jpg", "png"):
        raise ValueError("format must be jpg or png")
    if filter not in list("grizy"):
        raise ValueError("filter must be one of grizy")
    url = geturl(ra, dec, size=size, filters=filter, output_size=output_size, format=format)
    r = requests.get(url[0])
    im = Image.open(BytesIO(r.content))
    return im

def circle_points(RN, center, radius):
    radians = np.linspace(0, 2 * np.pi, RN)
    print("============")
    print(center[1])
    print(center[0])
    c = center[1] + radius * np.cos(radians)  # polar co-ordinates
    r = center[0] + radius * np.sin(radians)
    #print(np.array([c, r]).T)
    return np.array([c, r]).T


def cutImage(fullimage, xmin, xmax, ymin, ymax, otname, ra ,dec, rot):

    from datetime import datetime
    startTime = datetime.now()  

    newimg_data = []
    fullimage = fullimage
    print(fullimage)
    xmin = xmin
    xmax = xmax
    ymin = ymin
    ymax = ymax
    otname = otname
    ra = ra
    dec = dec
    rot = rot
    imgsize = 256
    #gim = getgrayim(ra, dec, size=imgsize, filter="g")
    gim = getcolorim(ra, dec, size=imgsize, filters="grizy")
    endTime = datetime.now()
    getUrlImgTime =(endTime - startTime).total_seconds()  
    #print(fullimage)
    #print(xmin)
    #print(xmax)
    #print(ymin)
    #print(ymax)

    # cutimage = "aa.png"
    cutimage = fullimage[0:-4] + "_" + otname + "_cut.png"
    cutimage2 = fullimage[0:-4] + "_" + otname + "_cut_2.png"
    cutimagefit = fullimage[0:-4] + "_" + otname + "_cut.fit"
    print(f' cutimage is {cutimage}')

    with pyfits.open(fullimage) as hdul:
        xprint_time
        #hdul.info()
        hdul.verify('fix')
        newimage = hdul[0]
        img_data = np.array(hdul[0].data)
        wcs = WCS(hdul[0].header)
        print(wcs)
        #print(img_data)
        #print(xmin)
        #print(xmax)
        #print(ymin)
        #print(ymax)
        # print(len(img_data))
        # print(img_data.shape)
        #newimg_data = img_data[xmin:xmax, ymin:ymax]
        newimg_data = img_data[ymin:ymax, xmin:xmax]
        # print(newimg_data)
        print("@@@@@@@@@@@@@@@@@")
       # print(len(newimg_data))
       # print(newimg_data.shape)
        
        pyfits.writeto(cutimage, newimg_data, header=None, overwrite=True)

        #================================================
        #to get the subimage in fit format
        #newimage.data = img_data[xmin:xmax, ymin:ymax]
        newimage.data = img_data[ymin:ymax, xmin:xmax]
        newimage.header = hdul[0].header
        phdr = newimage.header
        CCDSEC="[%d:%d,%d:%d]" % (xmin, xmax, ymin,ymax)
        BIASSEC="[%d:%d,%d:%d]" % (1,xmax-xmin+1, 1, ymax-ymin+1)
        TRIMSEC=CCDSEC
        #print(TRIMSEC)
        WCSDIM = 2
        LTM1_1 = 1
        LTM2_2 = 1
        WAT0_001='system=image'
        WAT1_001='wtype=tan axtype=ra'
        WAT2_001='wtype=tan axtype=dec'
        LTV1=-(xmin-1)
        LTV2=-(ymin-1)
        #print(LTV1)
        #print(LTV2)
        CRPIX1=phdr['CRPIX1']+LTV1
        CRPIX2=phdr['CRPIX2']+LTV2
        #print( CRPIX1)
        #print(CRPIX2)
        phdr.set('WCSDIM',WCSDIM)
        phdr.set('LTM1_1',LTM1_1)
        phdr.set('LTM2_2',LTM2_2)
        phdr.set('WAT0_001',WAT0_001)
        phdr.set('WAT1_001',WAT1_001)
        phdr.set('WAT2_001',WAT2_001)
        phdr.set('CRPIX1',CRPIX1)
        phdr.set('CRPIX2',CRPIX2)
        phdr.set('CCDSEC',CCDSEC)
        phdr.set('BIASSEC',BIASSEC)
        phdr.set('LTV1',LTV1)
        phdr.set('LTV2',LTV2)
        #phdr.set('TRIMSEC',TRIMSEC)
        #phdr.set('IMWHOLE',fullimage)
        phdr.set('EXTEND',False)
        #print(phdr)
        newimage.writeto(cutimagefit,overwrite=True)

        endTime = datetime.now()
        saveNewImgTime =(endTime - startTime).total_seconds()  

        #============================
        xprint_time
        interval = ZScaleInterval()
        vmin, vmax = interval.get_limits(newimg_data)
        #print(vmin)
        #print(vmax)
        points = circle_points(200, [50, 50], 15)[:-1]
        # plt.subplot(1, 1, 1)
        #fig = plt.gcf()
        #plt.axis('off')
        #fig.set_size_inches(7.0 / 3, 7.0 / 3)
        #plt.gca().xaxis.set_major_locator(plt.NullLocator())
        #plt.gca().yaxis.set_major_locator(plt.NullLocator())
        #plt.subplots_adjust(top=1, bottom=0, left=0, right=1, hspace=0, wspace=0)
        #plt.margins(0, 0)
        print("to show the wcs")
        plt.imshow(newimg_data, cmap=plt.cm.gray, vmin=vmin, vmax=vmax, origin='lower', aspect='equal')
        #print("to show the circle")
        plt.plot(points[:, 0], points[:, 1], color='white', linestyle="dashed", lw=1)
        #plt.plot(points[:, 0], points[:, 1], color='white', linestyle="dashed", lw=1)
        plt.savefig(fname=cutimage2, format='png', dpi=500)
        xprint_time
        #plt.subplot(211,projection=wcs)
        plt.subplot(211)
        print("to show the imshow")
        #im1 = plt.imshow(newimg_data, cmap=plt.cm.gray, vmin=vmin, vmax=vmax, origin='lower', aspect='equal')
        # rot = "80"   #定义了旋转角度，目前CGFT主焦点是旋转量，先上下镜像，再转转
        #rot = calRotZero(fullimage)
        #rot = rot
        im1 = pil.open(cutimage2)
        im2 = im1.transpose(pil.FLIP_TOP_BOTTOM)
        #im3 = im2.rotate(360 - eval(str(rot)))
        im3 = im2.rotate(0 - eval(str(rot)))
        plt.imshow(im3)
        #print("to show the circle")
        #plt.text(50, 50, otname, fontsize=20, color="white", verticalalignment="bottom", horizontalalignment="center")
        plt.title(f"OT name is {otname}",fontsize=10)
        plt.subplot(212)
        points = circle_points(100, [128, 128], 20)[:-1]
        plt.plot(points[:, 0], points[:, 1], color='white', linestyle="dashed", lw=1)
        print("to show the ps1 png")
        #plt.imshow(im1)
        plt.imshow(gim)
        print("to save the image")
        plt.savefig(fname=cutimage, format='png', dpi=500)
        xprint_time
        #plt.show()

        endTime = datetime.now()
        totalTime =(endTime - startTime).total_seconds()  

        plotImgTime = totalTime - saveNewImgTime
        saveNewImgTime = saveNewImgTime - getUrlImgTime

        print("cutImage total %.1f seconds, getUrlImgTime %.1f, save2NewImgTime %.1f, plotImgTime %.1f" % (totalTime, getUrlImgTime, saveNewImgTime, plotImgTime)) 


if __name__ == "__main__":
    #img = "G171013_C05906_0696.fit"
    # img = "data/sb22113056_g_1670072664_0001_1.fit"
    # xmin = 1421.96
    # xmax = 1521.96
    # ymin = 2568.18
    # ymax = 2668.18
    # otname = "aa"
    # ra = 30.40003
    # dec = 50.90454
    # rot = 30  #rotation
    img = str(sys.argv[1])
    xmin = str(sys.argv[2])
    xmax = str(sys.argv[3])
    ymin = str(sys.argv[4])
    ymax = str(sys.argv[5])
    otname = str(sys.argv[6])
    ra = str(sys.argv[7])
    dec = str(sys.argv[8])
    rot = str(sys.argv[9])
    print(img)
    print(xmin)
    print(xmax)
    print(ymin)
    print(ymax)
    cutImage(fullimage=img, xmin=int(float(xmin)),
             xmax=int(float(xmax)), ymin=int(float(ymin)),
             ymax=int(float(ymax)), otname=otname,
             ra=float(ra), dec=float(dec), rot=float(rot))

